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Abstract 

>' 

\^ • We present a formalism that leads very naturally to a hierarchi- 

cal description of the different contrast structures in images, provid- 
ing precise definitions of sharp edges and other texture components. 

C3 ' Within this formalism, we achieve a decomposition of pixels of the 

image in sets, the fractal components of the image, such that each set 

only contains points characterized by a fixed stregth of the singular- 

C^ ' ity of the contrast gradient in its neighborhood. A crucial role in this 

description of images is played by the behavior of contrast differences 
under changes in scale. Contrary to naive scaling ideas where the im- 
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pj ' age is thought to have uniform transformation properties ( Field 1987 ) 



O ■ each of these fractal components has its own transformation law and 

scaling exponents. A conjecture on their biological relevance is also 



given. 
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1 Introduction 

It has been frequently remarked that natural images contain redundant in- 
formation. Two nearby points in a given image are very likely to have very 
similar values of light intensity, except when an edge lies between them. In 
this case, the sharp change in luminosity also represents a loss of predictabil- 
ity, because probably the two points belong to two different objects, with 
almost independent luminosities. Since they cannot be predicted easily, it is 
frequently said that edges are the most informative structures in the image 
and that they represent their independent features ( |Bell fc Sejnowski||1997|) . 



Anyone's subjective experience is consistent with this: most people would 
start drawing a scene by tracing the lines of the obvious contours in the 
scene. Even more, the picture can be understood by anyone observing that 
sketch. Only afterwards, one starts adding texture to the light-flat areas in 
the drawing. These textures represent the less informative data about the 
position and intensity of the light sources. It seems that different structures 
can be distinguished in every image, following a hierarchy defined according 
to their information content, from the sharpest edges to the softest textures. 
Image processing techniques also emphasize the fact that edges are the 
most prominent structures in a scene. A large variety of methods for edge 
detection have been proposed, e.g. high pass filtering, zero-crossing, wavelet 
projections, ... ([Gonzalez fc Woods| |1992| ). But once the edges have been 



extracted from the image, these methods do not address the problem of 
how to describe and obtain the other structures. In order to detect these 
other sets and to extract them in the context of an organized hierarchy, 
it is necessary to provide a more general mathematical framework and, in 
particular, to define the variables suitable to perform this classification. 

In this paper we present a more rigorous approach to these issues, propos- 
ing precise definitions of edges and other texture components of images. As 
we shall see, our formalism leads very naturally to a hierarchical description 
of the different structures in the images, achieving a decomposition of its 
pixels in sets, or components, ranging from the most to the less informative 
structures. 

Another motivation to understand the role played by edges in natural 
images is related to the development of the early visual system. The brain 
has to solve the problem of how to represent and transmit efficiently the 
information it receives from the outside world. As it has been frequently 



pointed out ( [Barlow| PIMT] , |Llnsker| |T^ , \Mic^ [T^^ |van Haterenj |T^U^ , 
cells located in the first stages of the visual pathway do this by trying to 
take advantage of statistical regularities in the image ensemble. The idea is 
that since natural images are highly redundant over space, they have to be 
represented in a compact and non-redundant way before being transmitted 
to inner areas in the brain. In order to do this, it is first necessary to 
know the regularities of the images, which then could be used to build more 
appropriate and efficient internal representations of the environment. 

The statistical properties of natural images have been studied for many 
years, but the effort was mainly restricted to the second order statistics. 
The emphasis was put on the two-point correlation first because of the sim- 
plicity of the analysis, and secondly because its Fourier transform (i.e. the 
power spectrum) has a power law behavior that refiects the existence of scale 
invariant properties in natural images. This implies that the second order 
statistics does capture some of the regularities. However most of the under- 
lying statistical structure remains unveiled. One can easily convince oneself 
that most of the correlations in the image are still present after whitening is 
performed: the contours of the objects in the whitened image are still easily 
recognizable ( [l''ield| p.987| ). The conclusion is that non-gaussian statistical 
properties of edges are important to describe natural images, a fact that has 
also been recognized by ( [Ruderman fc Bialek] |1994| [Ruderman] |1994D . How- 
ever, no systematic studies of the statistical properties of edges had been 
done until very recently ( |l'uriel, Mato, Farga fc INadai 199^ ). 

A complete characterization of the regularities of natural images is im- 
possible, and probably useless. Some intuition is necessary to guide the 
search for regularities. The aim of this work is to exploit the basic idea that 
changes in contrast are relevant because they occur at the most informative 
pixels of the scene. Hence we focus our attention in the statistics of changes 
in luminosity. Since these changes are graded, the task could seem rather 
difficult and even hopeless, but the solution to this problem turns out to be 
quite simple and elegant. A single parameter is enough to describe a full 
hierarchy of differences in contrast and the hierarchy can be explicitly con- 
structed. Moreover, it is intimately related to the scale invariance properties 
of the image, which are certainly more elaborated than those that have been 
used until now but are still simple enough and have plenty of structure to 
be detected by a network adapting to the statistics of the stimuh. 

To understand these scale invariances we need first the basic concepts 



of multiscaling and inultifractality ( b'alconer 1990 ), together with a simple 



explanation of their meaning for images. While in naive scaling there is a 
single transformation law throughout the image (which is then said to be 
a fractal), in the more general case there is no invariance under a global 
scale transformation but the space can be decomposed in sets with a dif- 
ferent transformation law for each of them (the image is then said to be a 
mult if r act al). 

To make this concept more precise, let us consider a simple situation 
where a set of images present a single type of changes in contrast: the lu- 
minosity suddenly changes from one (roughly) constant value to another. 
In this case, a natural description of the properties of changes in contrast 
would be to count the number of such changes that appear along a segment 
of size r, oriented in a given direction. The statistical distribution of the 
position and size of these jumps in contrast would give relevant information 
about an ensemble of such images. However, in real images this counting is 
not enough to describe well the properties of contrast changes. This is be- 
cause these changes do not follow such a simple pattern: some of them can 
be sharp (the notion of sharp changes has to be given), but there will also 
be all types of softer textures that would be lost by just counting the most 
noticeable changes. It is then necessary to consider all of them as a whole. 
A natural way to deal with contrast changes is by defining a quantity that 
accumulates all of them, whatever their strength, contained inside a scale r, 
that is: 

e,(f) = — / dx[ / dx' I VC(f' ) I . (1) 

Hereafter, the contrast C{x) is taken as C{x) = I{x) — (/), where I{x) is 
the field of luminosities and (/) its average value across the ensemble. The 
bidimensional integral [] on the right hand side, defined on the set of pixels 
contained in a square of linear size r, is a measure of that square^ . It is 
divided by the factor r^, which is the usual Lebesgue measure (which we 
will denote A) of a square of linear size r. The quantity ej,(x) can then be 
regarded as the ratio between these two measures. More generally, we can 

^ The variables x'l , x'2 are the components of the vector x' and 1 1 denote the modulus 
of a vector. 

^ Any measure of a set has the property of being additive, that is: if one splits the set 
in pieces, the measure of the set is equivalent to the sum of the measures of all its pieces. 



define the measure /i of a subset A of image pixels as 

/i(A) = f dx' \ VC(x') I , (2) 



lA 

(/^ dx' means bidimensional integration over the set A) what we will call 
the Edge Measure (EM) of the subset A. It is also possible to generalize the 
definition of e^ for any subset A: 

eA is a density, a quantity that compares how much the distribution of |VC| 
deviates from being homogeneously distributed on A. We shall call it the 
Edge Content (EC) of the set A. Denoting Br{x) the "balF'Q of radius r 
centered around x, it is clear that er{x) = eB^(s), so this definition generalizes 
the previous one. We shall refer to er{x) as the Edge Content of x at the 
scale r. Notice that edr{x) is the /x-density at x, |VC|(a^). Thus, eq. (|[) can 
be symbolically expressed as 

d^{x) = \VC\{x)dx . (4) 

The main point is that contrast changes are distributed over the image in 
such a way that the Edge Content has rather large contributions even from 
pixels that are very close together. As a consequence, the measure presents 
a rather irregular behavior. But it is precisely in its irregularity where the 
information about the contrast lies. Let us assume for instance that, as r 
becomes very small, the measure behaves at every point x as 

lx{Br{x)) = /i,(x) = a(f) r^("")+'^ , (5) 

where rf = 2 is the dimension of the images. For convience a shift in the 
definition of the exponent of iXr^x) has been introduced. The Edge Content 
then verifies : 

tr{x) = a{x) r''^^^ , (6) 

and this choice of h{x) removes trivial dependencies on r. The exponent 
h{x) is just the singularity exponent of |VC|(x). That is, h < indicates 

^ "Ball" in a mathematical sense: it does not need to be true circles, but could also be 
squares or diamonds. In fact, along this paper the "balls" will be taken as squares. 



a divergence of |VC|(a;) to infinity, while h > indicates finiteness and 
continuous behavior. The greater the exponent the smoother is the density 
around that point. In general, we will talk about "singularity" even for the 
positive h cases. 

Thus, the meaning of eq. (^ is that all the points are singular (in this 
wide sense), and that the singularity exponent is not uniform. This irregular 
behavior would allow us to classify the pixels in a given image: the set of 
pixels with a singularity exponent contained in the interval [h — Ah, h + Ah] 
(where Ah is a small positive number) define a class Fh. These classes are 
the fractal components of the image. The smaller the exponent the more 
singular is the class, and the most singular component is the one with the 
smallest value of h. A measure fi verifying eq. (^ in every point x is called 
a multifractal measure. 

We will refer to these components of the images as fractal sets because it 
will be seen that they are of a very irregular nature. One way to characterize 
the odd arrangement of the points in fractal sets is by just counting the 
number of points contained inside a given ball of radius r; we will denote by 
Nr{h, Ah) this number for the set Fh- As r ^ it is verified that 

Nr{h,Ah) ~ r^e^). (7) 

This new exponent, D{h), quantifies the size of the set of pixels with sin- 
gularity h as the image is covered with small balls of radius r. It is the 
fractal dimension of the associated fractal component Fh, and the function 
D{h) is called the dimension spectrum (or the singularity spectrum) of the 
multifractal. It is also worth estimating the probability for a ball of radius 
r to contain a pixel belonging to a given fractal component. Its behavior for 
small scales r is: 

Prob{h,Ah) ~ r'^-^^^^ , (8) 

where d = 2 as before. Notice that, since r is small, this quantity increases 
as D{h) approaches d^ 

The multifractal behavior also implies a multiscaling effect: as every 
fractal component possesses its own fractal dimension, each of them changes 
differently under changes in the scale. Thus, although every component 

■* It can be proved that any set contained in a d-dimensional real space has fractal 
dimension D smaller or equal to d. 



has no definite scale and is invariant under the scale transformations, scale 
invariance is broken for a given image: the points in the image transform 
with different exponents, according to the fractal component to which they 
belong. This leads to the presence of intermittency effects: if the image 
were enlarged, its statistical properties would change. Self-similarity is lost. 
However, scale invariance is restored in the averages over an ensemble of such 
images. 

At this point one may wonder if a multifractal behavior is actually real- 
ized in natural images. One way to check its validity is to verify that eq. (^) 
is indeed correct. This in turn leads us to the technical problem of how to 
obtain in practice the exponent h{x) at a given pixel x of the image. This is 
a local analysis where the singularity is found by looking at the neighbour- 
hood of the point x. It has the advantage that it not only provides us with 
a tool to check the multifractal character of the measure but it also gives a 
way to decompose each image in its fractal components. These issues are 
addressed in the next Section. The multifractality of the measure can also 
be studied by means of statistical methods. In this case one looks for the 
effect of the singulaties on the moments of marginal distributions of, e.g., 
the Edge Content defined in eq. ([1|). This is done in Section]^, where the im- 
portant notions of self-similarity, extended self-similarity and multiplicative 
processes are explained. Once the relevant mathematical tools are given, 
these are used in Section ^ to analyze natural images. First in subsection 
|4.2| we check that the measure introduced in eq. (^ is indeed multifractal. 
Examples of fractal components of natural images are also presented in the 
same section. After this we study the consistency between the singularity 



analysis and the statistical analysis of the measure (subsection |4.3| ). The re- 
sults of the numerical study of self-similarity properties of the contrast C{x) 
and of the measure density |VC|(a;) are presented in subsection ^^ . The 
discussion of the results and perspectives for future work are given in the 
last section. Some technical aspects are included in the Appendix. 

2 Singularity analysis: The wavelet transform 

In this section we present the appropriate techniques to deal with mul- 
tifractal measures, and especially with those adapted to natural images. 
To check if /i defines a multifractal measure and to obtain the local ex- 
ponent h{x), we need to analyze the dependence of fj,{Br{x)) on the scale 



r. Unfortunately, direct logarithmic regression performed on eq. (^ yields 
rather coarse results when it is applied on discretized images, allowing only 
the detection of very sharp, well isolated features. It is then necessary to 
use more sophisticated methods. A convenient technique for this purpose 
( IMallat fc Huang||1992| , |Arneodo, Argoul, Bacry, Elezgaray fc Muzy| |1995| , |Ar^ 



neodo|p"996D is the wavelet transform (see e.g. ( |Daubechies| |1992| )). Roughly 



speaking, the wavelet transform provides a way to interpolate the behavior of 
the measure /i over balls of radius r when r is a non-integer number of pixels. 
Besides, if the image is a discretization of a continuous signal, the wavelet 
analysis is the appropriate technique to retrieve the correct exponents, as it 
is explained later on. 

The wavelet transform of a measure /i is a function that is not localized, 
neither in space nor in frequency. Its arguments are then the site x and the 
scale r. It also involves an appropriate, smooth interpolating function \I'(x). 
More precisely, the wavelet transform of /x is: 

TLd^ix) = /c//i(f')^r(x-f) = -^®^^(f) , (9) 

J ax 

where \l/r(a^) = ^^(-)- The important point here is that for the case of a 
positive measure it can be proved that: 

ndfi{xo) ~ r'^(^^o) (10) 

if and only if fi verifies eq.(^ {h being exactly the same) and \1/ decreases 
fast enough (see ( [Arneodo| |199(j| ) and references therein). In that way, this 



property allows us to extract the singularities directly from T^dfx. It is 
very important to remark that \l/ can be a positive function, what makes an 
essential difference with respect to the analysis of multiaffine functions (see 
section ^.p| 

In general the ideal, continuous measure dfi(x) is unknown. In these 
cases the data is given by a discretized sampling, which in our case is a 
collection of pixels. One can reasonably argue that the pixel is the result 
of the convolution of an ideal signal with a compact support function Xa{x) 
(describing, e.g., the resolution ~ a of the optical device). In this case, under 
convenient hypothesis, the wavelet transform yields the correct projection of 



^ This means that the exponents can be obtained with a non- admissible wavelet ^, 
i.e. a function with a non-zero mean 



the ideal signal. Let us call X(^ni,n2) — ^n the positions of the pixels in 
the discretized sample and C^ the values of the discretized version of the 
contrast; then 

Ch = T^C{Xn) =C(» XaiXn) , (H) 

Taking the sample {x,?} in such a way that the centers of adjacent pixels are 
at a distance larger than a, it is verified that: 



and so the discretized version fj,^'^^^> of the ideal measure /i, can be considered 
as a wavelet coefficient of the continuous fi with the wavelet x ^^ the scale 
a, namely 

dfif'^ ~ dx |VC| (g) Xaixft) = dx T^dfiixfi) ■ (13) 

Taking into account that the convolution is associative, the discretized wavelet 
analysis of d/x^ *^ with a wavelet \l/ can be simply expressed as: 

nd^i^^'^^xn) ~ n^^^ dfi{xn) (14) 

that is, analyzing the discretized measure /i'-'^*'^^ with a wavelet \l/ is equivalent 
to analyze the continuous measure /x with a wavelet \I/®X-- But if the scale r 
is large enough compared to the photoreceptor extent a, X-{^) ~ ^{x), which 
does not depend on r. Thus, if the internal size a is small enough to allow 
fine detection, we can recover the correct exponent at the point x^j by just 
analyzing the discretized measure. Conversely, the extent a imposes a lower 
cut-off in the details that can be observed in any image constructed by pixels: 
as it is just the size of the pixel, the estimation of the singularity exponent at 
a given point requires to perform the analysis over scales containing several 
pixels. 

This analysis is similar to that performed in ( pVIallat fc Zhong| |1991| , |Mal- 



lat fc HuangI p!992D : any image is analyzed by means of its scaling properties 
under wavelet projection. There are however important differences between 
the study performed by these authors and the present work. The first dif- 
ference concerns the basic scalar field to be analyzed: Mallat and Zhong 
considered C{x) instead of |VC|(a;), which we prefer because, as it will be 

9 



shown later on, it is the /i-density of a multifractal measure fi (see Section ^. 
Another difference is methodologicah in this paper we intend to classify the 
points according its singularity exponent to form every fractal component, 
while the previous papers were devoted to obtain only the wavelet transform 
modulus maxima, a concept related but different from what we will call the 
most singular component, that is, only one (although the most important) of 
the fractal components. The third difference concerns motivation: we want 
to classify every point with respect to a hierarchical scheme, while these 
authors search a wavelet-based optimal codification algorithm, in the sense 
of providing an (almost) perfect reconstruction algorithm from the wavelet 
maxima. 



3 Statistical self- similarity of natural images 

In this section we introduce the basic concepts of the statistical approach, 
defining at the same time new variables which are closer to the wavelet 
analysis. 

3.1 SS and ESS 

From the statistical point of view, Hr is a random variable, in the sense 
that given a random point x in an arbitrary image X the value of [Ji^{x) 
cannot be deterministically predicted. 

By definition, the multifractal structure is a multiscaling effect. This mul- 
tiscaling should be somehow apparent from the dependence of the probability 
density of /i,. on the scale parameter r. To determine the whole probability 
density function of Hr requires a large dataset, especially for the rare events; 
besides its dependence on r could mix multiscaling effects in a complicated 
way. On the contrary, the p-moments of ^r (that is, the expectation values 
of /ir raised to the p-th power) suffice to characterize completely the prob- 
ability density,^ and given that they are averages of powers of /i^ they are 
likely to depend on r in a simple way, somehow related to the behavior of 
the measure in eq. (^). 

It is more convenient to deal with the p-moments of normalized versions 
of Hr such as the Edge Content e,. = ^/i^ or the wavelet projections \T^dn\. 

^ Provided they do not diverge too fast with p. 
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They are normalized in such a way that the first order moments, (e^) and 
(|T^(i/i|), do not depend on r: the trivial factor r^'^ is hence removed from 
the random variable. 

As we will see later on, the moments (e^) and {\T^dfi\P) exhibit the re- 
markable scaling properties of Self-Similarity (SS): 

{eP) = a^r^v (15) 

(analogously for |T^(i/i| with exponents r*) and of Extended Self-Similarity 
(ESS) (referred to the moment of order 2) (|Benzi, Ciliberto, Baudet, Chavar^ 
ria fc 'Iripiccioii^ |1993| , [Benzi, Ciliberto, 'Iripiccione, Baudet, Massaioli fc 



Succi| p!993|, Penzi, Ciliberto fc Chavarria) |1995 



{el) = A{p,2) {el) 



pip,'i) 



(16) 



(and similarly for \T^dfi\ with exponents p*(p, 2)). ESS has been referred 
to the moment of order 2, but eq ( [T^ ) trivially implies that any moment 
can be expressed as a power of the moment of order q with an exponent 

In this way, we have the sets of SS exponents Tp and ESS exponents p{p, 2) 
for er and similarly the exponents r* and p*(p, 2) for T^dfi. Notice that if 
SS holds also does ESS, in which case there is a simple relation between Tp 
and p(p, 2): 



Tp = ^2P(P,2) 



(17) 



The actual dependence of Tp and p{p, 2) on p determines the fractal struc- 
ture of the system. A trivial dependence Tp oc p would reveal a monofractal 
structure. A different dependence on p indicates the existence of a more 
complicated (multifractal) geometrical structure. 



3.2 The multiplicative process 

A random variable subtending an area of linear size r and which possesses 
SS (eq. (|15D) and ESS (eq. (p!6D) can be described statistically by means 
of a multiplicative process ( [Benzi, Biferale, Crisanti, Paladin, Vergassola fc| 
Vulpiani] |1993| , [Novikov||1994| ). The statistical formulation is rather simple: 
given two different scales r and L such that r < L, there is a simple stochastic 
relation between the Edge Contents at these two scales : 
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ttr-L CL 



(18) 



where a^L is a random variable, independent of the Edge Content, which 
describes how the change between these two scales takes place. The symbol 
"=" indicates that both sides of the equation are distributed in the same 
way, but this does not necessarily imply the validity of the relation at every 
pixel of a given image. The random variables a^L between all possible pairs 
of scales r and L are said to define a multiplicative process. 

Let us understand in a more intuitive way the meaning of this process. 
It is rather clear that given the distribution of e^ at a fixed large scale L, to 
compute the distribution of the Edge Content e^ it is enough to know the 
distribution of a^L- But let us now introduce an intermediate scale r'. One 
could also obtain first e^i using ar'L and then e^ by means of a^r'- Thus, the 
multiplicative process must verify the cascade relation: 



«rL 



Oirr' Oir'L i 



(19) 



that tells us that the variables UrL are infinitely divisible. Infinitely divisi- 
ble stochastic variables are then completely characterized by their behavior 
under infinitely small changes in the scale, that is, by ar^r+dr- The simplest 
example of an infinitely divisible random variable UrL is the one having a 
binomial infinitesimal distribution, namely: 



(^r,r-\-dr 



1-A^ 

r 

/5(1-A? 



with probability 1 — [d — D^o] 
with probability [d - D^o] — 



(20) 



The form of this process is determined by its compatibility with eq. ([19|); see 
for instance ( |She fc Levequej |1994| , ^he fc Waymire| |1995| , [Castainj |1996| ) for 
a detailed study. Moreover, there are only two free parameters: it is verified 
that d — DoD = A/(l — (3). They have also to satisfy the following bounds: 
< /? < 1 and < A < 1. 

The infinitesimal process ( PO] ) can be easily interpreted in terms of the 
Edge Content and the measure density: 
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• Probability = 1 — [c/ — Doo] — : this is tlie most likely situation. In 
this case the Edge Content (and consequently the measure) changes 
smoothly under small changes in scale: Or^r+dr is close to one and its 
variation is proportional to dlnr. 

• Probability = [d — Doc] — : in this unlikely case, under infinitesimal 
changes in scale the Edge Content undergoes finite variations, which is 
refiected in the fact that ar,r+dr deviates from one by a factor jS. The 
/i-density | VC| must then be divergent somewhere along the boundary 
of the ball of radius r, and the parameter /3 is a measure of how sharp 
this divergence is. 

The process for non-infinitesimal changes in scale can be easily derived 
from the infinitesimal process presented in eq (^). In this case the ran- 
dom variable a^z, follows a Log-Poisson process ( ^he fc Waymire||1995|) ; its 
probability distribution parLi^ri) has the form: 



ParLyO^rL) 



^0 ^^- L rj \LJ 



If a multiplicative process is realised the Edge Content has the SS and the 
ESS properties, eqs. (^) and (|I6]). Thus, knowledge of the ESS exponents 



p{p, 2) and of T2 is enough to compute Tp using eq. (0). For the Log-Poisson 
process the ESS exponents are: 

which depend only on (3. This process was proposed in ( ^he fc Leveque|p"994|) 
for turbulent fiows, so we will refer to it either as the She-Leveque (S-L) or 
as the Log-Poisson model. The SS exponents are: 

Tp = -Ap + {d-D^){l-pP) (23) 

which compared with eq. ( p!7D gives the following set of relations : 

~K! ^ ^ (24) 

(1-/3)2 i_^ 
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This means that the model has only two free parameters, that can be chosen 
to be T2 and /? or A and Doo, for instance. The geometrical interpretation 
of the last two is very interesting ( [luriel et al.||199^ ), and will be explained 



in the next subsection. For the time being let us notice that the maximum 
value of the EC, || e^ ||oo, also follows a power law with exponent —A, 

II f II = re r~^ (2F)) 

II "^r ||oo — "oo ' \^^) 

and thus the parameter A characterizes the most divergent behavior present 
in natural images. 

3.3 SS and multifractality 

Let us now analyze in more detail the geometrical meaning of the model and 
its relation with the singularities described in Section ^. For a multifractal 
measure /i with singularity spectrum D{h) there is an important relation 
between the SS exponents Tp and D{h). Let us denote by p/i(/i) the proba- 
bility density of the distribution of the singularity exponents h in the image. 
Then, partitioning the image pixels according to their values of /i, eq. (^) 
can be expressed, for r small enough, as: 

{el) ~ ldhp^{h){a^)^r^^+'-''^^^ . (26) 



The factor r'^-^'^^'-) comes from the probability of a randomly chosen pixel 
being in the fractal F^ (see eq. (|^)) and guarantees the correct normalization 
of {cx^) p (see e.g. ( [brisch] p.995| , [Arneodo] |1996| ) for details). For very small 



r the application of the saddle point method to eq. (^) yields: 

I p\ r„ min\ph + d — D(h)} /^„x 

This gives a very interesting relation between Tp and D{h): Tp is the Legendre 
transform of D{h). The important point is that D{h) is easily expressed in 
terms of Tp by means of another Legendre transform, namely: 

D{h) = mm{ph + d-Tp} . (28) 

If one has a model for the Xp's (the S-L model, eq. (^)) the dimension 
spectrum can be predicted: 
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m^D^-'-±^ 






(29) 



This is represented in Fig. ^j. There are natural cut-offs for the dimension 
spectrum; in particular, there cannot be any exponent h below the minimal 
value hoo = —A. Thus, A defines the most singular of the fractal com- 
ponents, which has fractal dimension D{h^) = D^. We could use these 
properties of the most singular fractal (its dimension and its associated ex- 
ponent) to characterize the whole dimension spectrum, which is in agreement 
with the arguments given in section p.2| .f| 

3.4 Multiaffinity 

Related to multifractality there exists the simpler but more unstable 
property of multiaffinity ([Benzi, Biferale, Crisanti, Paladin, Vergassola fc 



Vulpiani] |1993|). This is a characterization of chaotic, irregular scalar func- 



tions that somewhat generalizes the concepts of continuity and differentia- 
bility. 

First, let us give the concept of Holder exponent: A scalar function F{x) 
is said to be Holder of exponent hp{xQ) at a given point xq if for any point 
y close enough to x the following inequality holds: 

|F(y)-F(fo)| < A,\y-x^\^^^^^^ , (30) 

Aq being a constant depending on the point Xq. We define the Holder ex- 
ponent of F{x) at Xq as the maximum of the exponents hplxo) verifying 
eq. (|30|). Defining the Linear Increment (LI) of F{x) by a displacement vec- 
tor r as 6ffF{x) = |F(x-|-r) — F{x)\, if the function F has a Holder exponent 
hp{x) at X then: 

5r-F{x) ~ apixy^^"^^ (31) 

which is similar to eq. (P). In the same spirit, we say that a function is 
multiaffine if for every point x eq. (^) is verified. 

Multiaffinity is a more intuitive concept than multifractality, because the 
Holder exponents are good characterizations of Taylor-like local expansions 
of the function F{x). For instance, hpix) > means that the function is 

^ Remember that P can be expressed in terms of A and Doc'. (3=1 — A/{d — I?oo) 
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continuous at x, hp{x) > 1 implies that the function has a continuous first 
derivative and so on. It also works in the other sense: a function having 
exponent hp = —1 behaves as "bad" as a (5-function (See the Appendix for 
a brief tutorial about Holder exponents of simple functions). Even more: a 
function F has Holder exponent hp at a given point if and only if its first 
order derivatives have Holder exponent hp — 1 at the same point. 

Besides, multiaffinity implies multifractality in the following sense: given 
a measure /x, if its density ^ is a multiaffine function then /x is a multifractal 
measure having exactly the same exponents as ^ at every point .0 For natural 
images and taking F = |VC|, this fact can be expressed as: 

e, ~ T^dfi ~ Sr-\VC\ (32) 

in the sense that the three variables have the same statistical dependence on 
r. This relation can be also expressed in terms of the SS exponents. Let us 
denote by Tp the SS exponents of e^, by r* those of T^d^ and by r^*" those 
of 5^1 VC|. Then the statistical relation in eq. ( p2|) implies that: 

(33) 

We are also interested in the relation between Tp and the SS exponents t^ of 
the moments of Sf^. If C exhibits the multiaffine behavior shown in eq. (pT|) , 
recalling the meaning in terms of differentiability of the Holder exponents, 
we have that 5^|VC|(x) verifies: 

6AVC\{x) ~ a^^,^{x) r'o^'^-' (34) 

This can be statistically interpreted as: 

er ~ 6AVC\ ~ -6r-C (35) 

r 

which refiects the shift by -1 in the exponents of the derivative of C.0 In 
terms of the SS exponents, this implies that: 

^ Notice that there is a shift of size d in the definition of the exponents of fi (eq. (||)) 
which is not present in the definition of the Holder exponents of a multiaffine function 
(eq. (1^)). 

^ Let us note that the relation rcr ^ SpC is the analog of the Kolmogorov hypothesis 



of local similarity (Arncodo 1996) 
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r^ = r^'^+P = Tp + p (36) 

Let us remark that the singularity spectra associated to the variables 6^^, 
(5^1 VC I and T^dfi are the Legendre transforms of the SS exponents r^ , t"^^ 
and r*, respectively (a fact that can be derived in the same way as eq. (]28|)). 
There still remains the converse question: if multifractality implies mul- 
tiafiinity. Unfortunately, this is not true (see ( |Daubechies| |199^ )). In fact, 
the \F{jJ) — F(xo)| term in eq. ( ^Ol) suggests the existence of a pseudo- Taylor 
expansion, in the sense that: 

F{y) ^ Pn{y-Xo) + Ao|y-Xo|"^("^°) (37) 

where Pn{y — ^q) is a polynomial of degree n in the projection of the vector 
y — Xq along its direction. We should thus generalize the concept of Holder 
exponent: given a point xq, we will say that this point possesses a Holder 
exponent n < hplxo) < n + 1 {n integer) if there exists a polynomial Pn{x) 
of degree n and a constant Aq such that 

\F{y)-p,^{y-x^)\ < A,\y-xo\''^^^'^^ (38) 

and hp^xo) is the maximum of the exponents verifying such relation. 

Multifractal measures with no multiafiine densities (as defined in eq. 
(|30D) have been studied in (Pacry, Muzy fc Arneodo| |1993| , |Arneodo||1996|) . 



The problem can be explained by the masking presence of the regular, poly- 
nomial part which stands for global, non-localized effects in the signal field 
F{y). To remove this part, those authors propose a generalization of the Lin- 
ear Increment. Instead of Linear Increments of the function, wavelet trans- 
forms of the appropriate one-dimensional restriction are considered, taking 
the scalar field as a real-valued measure density. This wavelet transform of 
the function F{x) is defined as 

r^F{xo) = [ dsF{xQ + s-)mr{s), (39) 

J r 

where \1/ is a real, one-dimensional analyzing wavelet. Since the integral is 
one-dimensional the wavelet at the scale r is ^r('S) = -\l/(-). If the scale r 
is small enough, the integral is dominated by the first terms in the pseudo- 
Taylor expansion of F{x) around Xq, eq. (^7\): 
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nF{xo) ^ Jdspn{s)^r{s) + Ao J ds \s\''^'^'^'^^r{s) (40) 

If the analyzing wavelet is now required to vanish at least the first n in- 
teger moments, the singular part appears as the main contribution to the 
wavelet transform. (For later reference, let us say now that a wavelet with 
a vanishing zero-th order moment - i.e., its integral is null - is called an 
admisible wavelet). Hence, changing variables t = - and defining A = 
Ao Jdt |t|'^^(^o)^(t), one has: 

T|F(xo) ^ Ar'^^(^o) (41) 

where all the dependence on r stands in the prefactor r^Fi^o) ^ jj^ this way, 
the Holder exponent (in the sense of eq. (^)) can be obtained by a wavelet 
projection, provided the wavelet has zero moments up to a large enough 
order. Since hpixo) is non- integer, this procedure cannot cancel the singular 
part even if the order n of the largest vanishing moment of \1' is larger 
than hp{xo). On the other hand, the wavelet could not be able to detect 
singularity exponents hp > n + 1. This is because a polynomial of order 
n + 1 could still be present and dominate the contribution of the singular 
term. For this reason, it is important to determine the right order of the 
required zero moments of the wavelet. 

This is then the appropriate scheme to detect the Holder exponents of 
signals with a regular part. In cases in which the /i-density |VC| and/or its 
primitive C are not multiaffine, the concept of multiaffinity (eq. (|30|) ) should 
be applied to their wavelet projections. It is then necessary to determine 
how many zero moments the wavelet should have, as it provides information 
on the global aspects of the signal. This generalization is straightforward, 
replacing in all the cases S^^ by T^C and (5^|VC| by T^|VC|. For instance, 
if the wavelet transforms of C and |VC| verify eq. (|4T| ) for certain A and 
hpix), the same arguments that led us to eqs. (|35[) and ( p6D can be used to 
derive that : 

er ~ T||VC| ~ -T|C (42) 

r 

which is expressed in terms of the SS exponents Tp ^ (of T^C) and Tp ^'^ (of 

r||VC|) as: 
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Tp 



Tp^^ +P = Tp+p 



(43) 



Let us make a final remark: it is possible that the wavelet projections 
of a given scalar function F verifies the scaling in eq. (^) without verifying 
eq. (pSl). This can happen, for instance, for functions with such an irregualar 
behavior that the pseudo- Taylor expansion is not possible. In those cases, 
T^F can be used as a stronger generalization of the usual Linear Increment. 



As it will be shown in Section [4. 41) , this is the case for |VC| 



4 Numerical analysis of images 

4.1 Databases 

We have used two different image datasets: one from Daniel Ruderman 
(see ( [Rudernian] |1994| )), consisting in 45 B/W images of 256 x 256 pixels 
with up to 13 bits in luminosity depth (from now on this set will be referred 
to as the first dataset); and another containing 25 images taken from Hans 
van Hateren's database (see (|van Hateren fc van der Schaaf] p!998|) ), having 
1536 X 1024 pixels and 16 bits in luminosity depth (we will refer to these 
images as the second dataset). These 25 images are enumerated in table |I| 



imk00034.imc 


imkOOSOl.imc 


imk01164.imc 


imk02649.imc 


imk03807.imc 


imk00211.imc 


imk00808.imc 


imk01406.imc 


imk03134.imc 


imk03842.imc 


imk00478.imc 


imk00881.imc 


imk02035.imc 


imk03514.imc 


imk03940.imc 


imk00586.imc 


imk01017.imc 


imk02603.imc 


imk03536.imc 


imk04042.imc 


imk00605.imc 


imk01032.imc 


imk02626.imc 


imk03789.imc 


imk04069.imc 



Table 1: The 25 images from van Hateren's database used in this 
study. They can be observed and downloaded from the URL address 
http://hlab.phys.rug.nl/archive.html 



The first dataset has been used for the statistical analysis of the one- 
dimensional wavelet transform, eq. (^), while the second dataset was used 
for the more demanding task (from the statistics point of view) of com- 
puting the moments of the two-dimensional variable in eq. ([^). No matter 
the ensemble considered, however, the studied properties seem to behave in 
essentially the same way. 
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4.2 Multifractality of the measure and singularity anal- 
ysis 

We first address the issue of tlie multifractality of the measure /i defined in 
eq. @. Our numerical analysis gives an affirmative answer to this question: 
the wavelet transform T^dfx{xo) is well fitted by ~ r'^(^°^ in about 98% of 
the pixels Xq in each image of our two datasets. The analysis was done 
using several wavelets, although the most convenient family of one parameter 
functions was \E''''(x) = q , i^pw (7 > !)• These are positive (non- admissible) 
wavelets. The singularity exponent of the measure at a given pixel, h{x), was 
obtained from a regression on eq. (0). The value of the scale r was typically 
taken in the range (0.5, 2.0) in pixel units. Notice that for scales smaller 
than r ~ 0.5 the method would detect an artificial "edge" that actually 
corresponds to the pixel boundaries, what corresponds to h = —1. Similarly, 
as the scale becomes large the wavelet can only discriminate between very 
sharp edges and, for r of the order of the size of the image it would yield the 
value h = everywhere. 

To obtain the fractal components Fh, the pixels of a given image were 
classified according to the value of their exponents. The exponents lie in the 
interval (—0.5, 1.4): only a negligible fraction of pixels have exponents less 
than —0.5 ( about 0.8% of them ) and none with a singularity below — 1. 
It is also observed that there are no exponents larger than 1.4. This range 
of possible exponent values coincides exactly with that predicted in ([Turiel 



et alj 1998 ), and is in agreemente with the model explained in section p.3| . 



The most singular component is the one with h = —0.5, defined within 
a window of size Ah. The visualization of this set is very instructive (see 
Fig. 0). Clearly, the points appear to be associated with the contour of the 
objects present in the image. This is in agreement with our expectation that 
the behavior of the measure at these edges should be singular, and as much 
as possible. But there is something rather surprising: for a sharp, theta- 
like contrast, with independent values at both sides of the discontinuity, one 
would have expected an exponent —1, which appears as a Dirac 5-function 
in the density measure VC (see example number 4 in the Appendix). Our 
interpretation of the minimal value of h, hoc = —0.5, is that it refiects the 
existence of correlations among different fractal components of the images 
that smoothen the singularities (in particular those of the sharpest edges). 
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It is also interesting to observe the fractal components with exponents 
just above the most singular one. Fig. |^ shows the two next components, 
defined with A/z. = 0.05. A comparison with Fig. |1| shows that pixels in these 
components are close to pixels in the most singular one. It is also rather clear 
that the probability that a random ball of size r contains a pixel from the 
component F^ increases with /i, at least for small values of this exponent. 
Since this probability scales as shown in eq. (H), this behavior reveals an 
increase with h of the dimensionality of the fractal components. It is also 
observed (although not shown in the figures) that the dimensionality starts 
decreasing beyond h ^ 0.2. This is in agreement with the arguments of 
section ^]^: the singularity spectrum D{h) of the S-L model, which is shown 
in Fig. ^, refiects this behavior of the fractal components. 

The analysis done in this subsection justifies the arguments presented in 
Section 0. One cannot consider the sharpest edges of an image independently 
of the other texture structures (i.e. the other fractal components). Luminos- 
ity jumps from one roughly constant value to another (statistically indepen- 
dent) constant value are not typical in natural images. The sharpest edges 
appear correlated with smoother textures, an effect that tends to decrease 
the strentgh of their singularity exponent, reaching the value h ~ —0.5. 
Standard edge-detection tecniques ( Ponzalez fc Woods] p.992|) have not been 



designed to search for singularities of this strength. Good edge-detection 
algorithms should be specifically constructed to capture singularities close 

to /loo- 

4.3 Statistical analysis of the measure 

We computed the p-moments of three different variables related to the mea- 
sure: 

• Bidimensional Edge Content (cr): Following the definition given 
in eq. (p. This analysis was done using the 25 images of the second 
dataset. As the variable itself is positive, we computed directly the 
quantities (e^). 

• Two one-dimensional restrictions of the wavelet projection of 
the measuref^ (T^dfi): We considered horizontal (h) and vertical (v) 



10 



Since the calculations involving T^dfj, are highly computer-time consuming we did 
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restrictions of the measure d^, that is dfi'^ = dx\-^\ and d/j,"" = dy\-T^\. 
The appropriate wavelet transforms are denoted by T^dfj!" {I = h,v). 
As these coefficients need not to be positive (in general \1/ is not posi- 
tive), we computed the moments of their absolute values, (|T^(i/x'|P). 

It is remarkable that the bi-dimensional Edge Content and the two one- 
dimensional wavelet projections exhibit the scaling properties of SS (eq. |T5|) 
and of ESS (eq. 0) (the ESS tests for the three variables are presented in 
Fig. H and H). This was also observed in ([Turiel et al.||1998D , although for 



somewhat different variables: the one-dimensional horizontal and vertical 
restrictions of the local edge variance. Surprisingly, no matter the particular 
variable, the exponents Tp and p{p, 2) obtained are very similar (see Fig. ^ and 
those in ( luriel et al. 1998 )). Therefore, they should refer to very essential. 



robust aspects of images. 

The S-L model fits very well the experimental exponents (see Fig. ^. 
It seems then that the multifractal structure underlying the statistical de- 
scription could be described with the Log-Poisson process, which means that 
under infinitesimal changes in the scale there are only two possible different 
types of transformations. 

It is really remarkable that for our datasets A ^ 0.5 and Dqo ~ 1- This 
means that the most singular component is a collection of lines (the sharp 
edges present in the image) which can be characterized by their common 
exponent —A = —0.5. But this is precisely the smallest value we observed 
previously in our wavelet analysis. The statistical approach confirms the 
result obtained by a local singularity analysis. 

The non-linear dependence of p{p, 2) on p again indicates that the Edge 
Measure is multifractal. This has to be contrasted with the results by 
(IRuderman fc Bialek 1994 , [Ruderman 1994|) who find that the log-contrast 



distribution follows a simple classical scaling. |^ It also holds for any derived 
variable as the gradient (see ( [Rudernianl |1994]) for details). This particular 



scaling is related, in that context, with general properties of scale invari- 
ance of natural images: the averaging of the contrast does not keep track of 
the details (e.g., the edges completely contained into the averaging block) 

not consider the bi-dimensional wavelet transform. FFT cannot be safely used due to 
aliasing, which heavily changes the tails (rare events) of the distribution. 

^^ Other luminosity analysis performed by D. Ruderman provided some evidence of 
multiscaling behaviour (private communication). 
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and produces a statistically equivalent image. The Edge Content and the 
wavelet projections, being local averages of a gradient, are able to keep the 
information about the more complex underlying multifractal structure. 

The same authors have also designed an iterative procedure that tends 
to decompose the image into a gaussian-like piece (local fluctuations of the 
log-contrast) and an image-like piece (local variances). Such a decomposition 
would suggest the existence of a multiplicative process for the log-contrast 
itself where the random variable relating the image at two different levels of 
the hierarchy is given by the local fluctuations. A closer look to the iterative 
procedure indicates that this is not the case: independence between the two 
factors of the process is not guaranteed and also the process is not infinitely 
divisible. 

4.4 MultiafRnity properties of C and |VC| 

In this section we study the scaling properties for the Linear Increments 
of the contrast C{x) and of the measure density |VC|(x). Our aim is then 
to check if C and |VC| satisfy the scaling defined in eq. (0). Recalling 
the connection (the analog of eq. ( p6| ) for these two quantities) between 
local scaling (described by the singularity exponents) and statistical scaling 
(characterized by the SS exponents), one concludes that the validity of SS 
implies that multiafiinity holds. We have then done a statistical analysis of 
the Linear Increments of the contrast and the measure density. 

We first computed the p-moments of the Linear Increments of both vari- 
ables, looking for SS and ESS scahngs. Let us notice that there is no a 
priori argument to expect that (5^(x) shows SS. On the other hand, since 
the measure /i is multifractal (as it was shown in section [4.2|) it would be 



reasonable to expect that |VC|(a;) has a multiafiine behavior. However, this 
has to be verified explicitly because the multifractality of the measure does 
not guarantee the multiafiinity of its density, as it was discussed in section 



3.4 . In fact, our numerical analysis showed that neither C{x) nor |VC|(x) 
have SS. The same is true for ESS, as it is shown in Figs, ^a and 0, where 
the moments exhibit a rather erratic behavior. 

It seems that the failure of SS for C{x) and |VC|(a;) has to be explained 
in different terms. In the case of |VC|(x), this quantity is the measure den- 
sity and, as it was stated in eq. (0), its convolution with a non-admissible 
(i.e. positive) wavelet has the same singularities as the measure, eq. (El). 
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This implies that the lack of SS cannot be explained in terms of a polyno- 
mial part as described in eq. (pTl). In fact, a positive wavelet cannot alter 
the contribution from this polynomial while one observes that the data for 
(5^-|VC|(x) (Fig. 0) are different from those of T^dfi^ and T^dfi" (Fig. |). 

On the other hand it is plausible that the failure of SS of 6i^{x) can 
be attributed to the presence of a polynomial part in C{x). Differently to 
what happened with |VC|(x), when the contrast is convolved with a non- 
admissible wavelet, SS and ESS are still not present (that is, one observes 
a similar erratic behavior as in Fig. ^a). The same behavior is observed 
when C{x) is convolved with an admissible wavelet with a zero mean. The 
situation changes drastically when a wavelet with vanishing zero- and first- 
order moments is used. The result obtained using the second derivative of 
a gaussian is presented in Fig. ^d, which clearly shows ESS (it also exhibits 
SS). This can be understood recalling that the possible values of the Holder 
exponent of C{x) are contained in the interval (0.5, 2.4) (see eq. (|3^ ) and 
section ^^ ). A non-admissible wavelet truncates the detected singularities 
at the value /i = 0, an admissible one truncates them at h = 1. The third 
wavelet considered, with two vanishing moments, truncates the singularities 
only from h = 2 and then it covers almost the whole interval (In theory a 
wavelet with another vanishing moment would be more appropriate, however 
in practice it is numerically more unstable because the resolution of the 
wavelet is smaller). 

Finally, we also checked that the experimental values of Tp, Tp ^ and Tp ^^ 



satisfy the relation given in eq. (|43|). 

5 Discussion 

We have shown that in natural scenes pixels are arranged in non-overlapping 
classes, the fractal components of the image, characterized by the behavior of 
the contrast gradient under changes in the scale. In each of these classes this 
quantity exhibits a distinct power law behavior under scale transformations, 
what indicates that there is no a well-defined scale for the isolated compo- 
nents but that at the same time the scene is not scale invariant globally. 

This result implies that a given image is composed of many different 
spatial structures of the contrast associated to the fractal components, which 
can be arranged in a hierarchical way. In fact, how sharp or soft a change 
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in contrast is at a given point can be quantified in terms of tlie value of the 
scaling exponent at that site. The smaller the exponent the sharper the local 
change in contrast is. 

More precisely, we have dealt with two basic quantities: the contrast 
C{x) and the edge measure /i, which density is the modulus of the contrast 
gradient |VC|(x) (see eq. (|^)). Closely related to the measure, we can define 
the edge content of a region A of the image (eq. (^), which quantifies how 
much the contrast gradient deviates from being uniformly distributed inside 
the region A. This is a suitable definition to study the local behavior of the 
contrast gradient. As the size of A is reduced, some contributions of the 
contrast gradient to the edge content are left outside this region. Under an 
infinitesimal change in the scale very often these contributions are small, but 
sometimes it does happen that a small change in the scale produces a large 
change in the edge content. This is due to singularities in the contrast gra- 
dient. As a consequence of this, the edge content (and the measure) exhibits 
large spatial fluctuations. The edge content can have different singularities 
at neighbouring points, characterized by a power-law behavior with expo- 
nents h{x) which depend on the site, eq. (^. This means that the measure 
is multifractal, as it has been explicitly verified in this work. 

The set of pixels with the same value of h{x) defines a particular spatial 
structure of contrast: the fractal component F^. The smallest h{x) is the 
sharpest the contrast gradient at the pixel x is. Its minimal value was, for 
our dataset, /loo = —0.5. Another question refers to the fractal dimension 
Doo of this component. One finds D^o = 1 ( |Turiel et al.||1998D , and when this 



component is extracted from the image one checks that it corresponds to the 
boundaries of the objects in the scene. The associated contrast structure is 
made of the sharpest edges of the image. The other, less singular components 
represent softer contrast structures inside the objects. It is also observed that 
there are strong spatial correlations between fractal components with similar 
exponents. 

We have also shown how to extract the fractal components from the 
image. To do this one needs a suitable technique to compute the scale 
behavior at a given pixel. A wavelet projection ( |Mallat fc Huan^ 1992 , 



|Arneodo| |1996| ) seems to be the correct way to perform this analysis, since it 



provides a way of interpolating scale values which are not an integer number 
of pixels. Our approach differs from the one followed in those references in 
several aspects, but in particular in that we have been here interested in the 
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characterization of the whole set of fractal components. 

We remark the important fact that the multifractal can be explained by 
a model with only two parameters (e.g. (3 and A). This had been noticed 
before in terms of a variable that integrates the square of the derivative of 
the contrast along a segment of size r, both taken either in the horizontal 
or in the vertical direction (luriel et al. 1998) . This quantity admits the 



interpretation of a local linear edge variance. We have seen that the same 
model is valid for many other quantities, some of them more general and 
natural. 

The basic idea of the model is that the contrast gradient has singularities 
that under an infinitesimal change in the size of A produce a substantial 
change in the edge content inside that region. The simplest stochastic process 
to assign a value to this change is a binomial distribution, eq (^Of ): its more 
likely value is very small, but with small probability this change is given by 
the parameter /3, which for our dataset is about 0.5. The other parameter 
of the model. A, measures the strength of the singularity. The geometrical 
locus of these singularities is the most singular fractal component. For a 
finite scale transformation this process becomes Log-Poisson. The model 
appeared before in other problems ( |She fc Leveque|[1994| ). 



In a sense, the most singular component conveys a lot of information 
about the whole image. The two model parameters can be obtained by 
analyzing properties of this particular class of pixels: its dimension (L'oo) 
and the scale exponent of the edge content (A = — /loo)- In turn, these 
two parameters completely define the whole dimension spectrum in the Log- 
Poisson model. It is then very plausible that other components have a great 
deal of redundancy with respect to the most singular one. 

This hierarchical representation of spatial structure can be used to predict 
specific feature detectors in the visual pathways. We conjecture that the most 
singular component contains highly informative pixels of the images that 
are responsible for the epigenetic development leading to the adaptation of 
receptive fields. Learning of the statistical regularities ( Parlow| |1961| ) present 



in this portion of visual scenes (using for instance the algorithm in ( [Bell 



fc Sejnowski||1997D ), would give relevant predictions about the structure of 



receptive fields of cells in the early visual system. 
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A Characterization of the singularities of sim- 
ple functions 



It is instructive to compute the singularity exponents defined in eq. ([30| ) for 
several simple functions. 

1. Continuous functions: If a function f{x) is continuous at a given 
point X, the value of / at y should be rather close to /(x) provided yis 
close to x; then |/(^| < |/(f)| +e so \f{x) - /(y)| < 2|/(f)| +e = A. 
That is, / is Holder of exponent at x. 

2. Smooth functions with continuous, non-vanishing first deriva- 
tive: Using the Taylor expansion, /(x) — f{y) = V/(^o) ■ {x — y) with 
yo a point between x and y. Let us call A = max|V/(yo)| > 0, so 
|/(x) — f{y)\ < A\x — y\. That is, / is Holder of exponent 1 at x. 

Analogously, any function / having n — 1 vanishing derivatives at a 
point X and a non-vanishing continuous n-th derivative, is Holder of 
exponent n at x. 

3. ^-function: The 6{x) function is given by 

r X < 

e{x) = < 1 X > 

[ undefined x = 

Since for x 7^ the increments are zero (for small enough displacement) 
it is Holder of any order. For x = we make use of the property that 
/ is Holder of exponent /i at a given point if and only if its primitive 
F is Holder of exponent h + 1 at the same point. 

One possible primitive F{x) of ^(x) is the following: 



F(x) 



X < 
X X > 



which is a continuous function. The increments \F{0) — F{r)\ = 
\F{r)\ = F{r) are thus immediate, and obviously F{r) < \r\. More- 
over, it is also clear that the maximal h such that \F{0)—F{r)\ < A\r\'^ 
is h = 1. Thus, F has Holder exponent 1 at x = and therefore 
9{x) has Holder exponent at x = 0. 
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4. (5-function: The Dirac's (5-function is in fact a distribution. Anyway, 
it is possible to define the Holder exponent of distributions attending 
to the fact that any distribution can be expressed as a finite order 
derivative of a bounded function. So, if the n-th order primitive of a 
distribution is a Holder function of exponent /i at a given point, the 
distribution is Holder of exponent h — n dX the same point. 

Since the 5{x) distribution vanishes in x 7^ the only non-trivial point 
is a; = 0. Since the 5 is the derivative of 9{x) and this has Holder expo- 
nent at a; = 0, one concludes that the 5{x) has Holder exponent —1 
at a: = 
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Figure 1: An image from Ruderman's ensemble ( [Rudernian| [1994|) and its 
most singular component (taken as the set of points having a singularity 
exponent in the interval h = —0.5 ± 0.05). 
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Figure 2: Fractal components with exponents h = —0.4 ± 0.05 (above) and 
h = —0.3 ± 0.05 (below) for the same image represented in Fig. |I[ Notice 
the spatial correlations between components with similar singularities. 
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Figure 3: Test of ESS for the bidimensional edge content e,.. The graphs 
represent the logarithm of its moments of order 3, 5 and 7 versus the loga- 
rithm of the second order moment, for distances r = 4 to r = 64 pixels, that 
is, about an order of magnitude smaller than the image size. This upper 
bound is needed both to consider values of r small enough compared to the 
size of the images and to be able to use moments of high order as p = 7. 
The data corresponds to the second image ensemble. According to eq. ( ]TB| ) 
these graphs should be straight lines; the best linear fits are also represented. 
Except for lower cut-off effects at very small distances, the scaling law is well 
verified. 
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Figure 4: Test of ESS for a) T^dfi^ and b) T^dfi"". The graphs represent 
the logarithm of the moments of order 3, 5 and 7 for both variables versus 
the logarithm of the corresponding second order moment. The scales range 
from r = 4 to r = 64 pixels. The first dataset was used in this computation. 
The corresponding best linear fits are also represented. Except for lower and 
upper cut-off effects at small and large distances, the scaling law verifies well. 
Although the test is not shown here, SS also holds for these variables. 
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Figure 5: ESS exponents p(p, 2) for a) the bidimensional edge content e,. , 
b) the horizontal wavelet transform T^dfi^ and c) the vertical one, T^dfi'^. 
The convolution function \1' was taken to be a gaussian function. Each 
value of p{p, 2) was obtained by linear regression of the logarithm of the 
p-th moment versus the logarithm of the second moment, for r between 8 
and 32. The solid line represents the fit with the Log-Poisson process. The 
best fit is obtained with a) (3 = 0.52 ± 0.05 b) /5 = 0.55 ± 0.06 and c) 
P = 0.5 ± 0.07 The SS parameters T2 were also calculated; they turned out 
to be a) T2 = -0.26 ± 0.06 , b) rg = -0.25 ± 0.06 and c) rs = -0.25 ± 0.08. 
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Figure 6: Test of ESS for a) the horizontal hnear increment of the contrast, 
b^-fi and b) the horizontal wavelet coefficient jT^C'*! (\E' was taken as the 
second derivative of a gaussian); the first dataset was used and the range of 
scales was from r = 4 to r = 64. ESS is not present for b^-fi ■, but it is clearly 
observed for the wavelet transform (notice the presence of the cut-offs due to 
numerical effects). Note that the scale of the axes in a) and b) is different. 
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Figure 7: Test of ESS for a) horizontal linear increment of |VC|, (5^^|VC| 
and b) vertical linear increment of the same function, 5^^|VC|. The first 
dataset was used and the scales were taken ranging from r = 4 to r = 64 
pixels. Again, these variables do not have ESS. The corresponding wavelet 
projections are shown in Fig. ^ where one can see that ESS holds. 
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Figure 8: Singularity spectrum of the Log-Poisson model for [3 = 0.5 and 
T2 = —0.25. For these values D^ = 1. and A = 0.5 (see text) 
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